#include "FLOAT.h"
#include <stdint.h>

FLOAT F_mul_F(FLOAT a, FLOAT b) {
    int result_high, result_low;

    asm volatile (
        "imull %2;"
        "shrdl $16, %%edx, %%eax;"
        : "=a"(result_low), "=d"(result_high)
        : "r"(b), "a"(a)
        : "cc"
        );

    return result_low;
}

FLOAT F_div_F(FLOAT a, FLOAT b) {
	/* Dividing two 64-bit integers needs the support of another library
	 * `libgcc', other than newlib. It is a dirty work to port `libgcc'
	 * to NEMU. In fact, it is unnecessary to perform a "64/64" division
	 * here. A "64/32" division is enough.
	 *
	 * To perform a "64/32" division, you can use the x86 instruction
	 * `div' or `idiv' by inline assembly. We provide a template for you
	 * to prevent you from uncessary details.
	 *
	 *     asm volatile ("??? %2" : "=a"(???), "=d"(???) : "r"(???), "a"(???), "d"(???));
	 *
	 * If you want to use the template above, you should fill the "???"
	 * correctly. For more information, please read the i386 manual for
	 * division instructions, and search the Internet about "inline assembly".
	 * It is OK not to use the template above, but you should figure
	 * out another way to perform the division.
	 */

    int result;
    int remainder;

    asm volatile (
        "movl %2, %%eax;"
        "sall $16, %%eax;"
        "cltd;"
        "idivl %3;"
        : "=a"(result), "=d"(remainder)
        : "r"(a), "r"(b)
        : "cc"
        );

    return result;
}

FLOAT f2F(float a) {
	/* You should figure out how to convert `a' into FLOAT without
	 * introducing x87 floating point instructions. Else you can
	 * not run this code in NEMU before implementing x87 floating
	 * point instructions, which is contrary to our expectation.
	 *
	 * Hint: The bit representation of `a' is already on the
	 * stack. How do you retrieve it to another variable without
	 * performing arithmetic operations on it directly?
	 */

    uint32_t float_bits;

    char float_bytes[4];
    *(float*)float_bytes = a;
    float_bits = *(uint32_t*)float_bytes;

    uint32_t sign = (float_bits >> 31) & 1;
    int32_t exponent = ((float_bits >> 23) & 0xFF) - 127;
    uint32_t mantissa = float_bits & 0x7FFFFF;

    if (exponent == -127 && mantissa == 0) {
        return 0;
    }

    mantissa |= 0x800000;

    int32_t result;
    if (exponent >= 16) {
        result = mantissa << (exponent - 16);
    }
    else if (exponent >= -16) {
        result = mantissa >> (16 - exponent);
    }
    else {
        result = 0;
    }

    if (sign) {
        result = -result;
    }

    return result;
}

FLOAT Fabs(FLOAT a) {
    return a < 0 ? -a : a;
}

/* Functions below are already implemented */

FLOAT sqrt(FLOAT x) {
	FLOAT dt, t = int2F(2);

	do {
		dt = F_div_int((F_div_F(x, t) - t), 2);
		t += dt;
	} while(Fabs(dt) > f2F(1e-4));

	return t;
}

FLOAT pow(FLOAT x, FLOAT y) {
	/* we only compute x^0.333 */
	FLOAT t2, dt, t = int2F(2);

	do {
		t2 = F_mul_F(t, t);
		dt = (F_div_F(x, t2) - t) / 3;
		t += dt;
	} while(Fabs(dt) > f2F(1e-4));

	return t;
}

